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We propose a phase model to study cascade failure in power grids composed of generators and 
loads. If the power demand is below a critical value, the model system of power grids maintains 
the standard frequency by feedback control. On the other hand, if the power demand exceeds the 
critical value, an electric failure occurs via step out (loss of synchronization) or voltage collapse. The 
two failures are incorporated as two removal rules of generator nodes and load nodes. We perform 
direct numerical simulation of the phase model on a scale-free network and compare the results with 
a mean-field approximation. 
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I. INTRODUCTION 



The dynamics in complex network has been intensively studied, since the proposal of the small world network 
and the scale-free network. 0, Q In the scale- free network, the degree distribution obeys a power law. Barabasi and 
Albert showed that a power grid in western United States was a typical example of scale-free networks. @ However, 
there is a controversy on the degree distribution of power grids. Some authors reported power-law distributions, 
whereas others reported exponential distributions. [H, 0] The distribution might be different for different countries and 
regions. In some regions, a two-dimensional planar network might be a good approximation. Furthermore, there is a 
hierarchical structure from the load dispatching center to power stations via many transformer substations in electric 
power systems. 

A serious issue in power grids is a large-scale failure of power supply, called a blackout. || A large-scale blackout or 
a cascade failure in a power grid network can have disastrous consequences in modern society. From a viewpoint of 
statistical physics, several simple models have been proposed to study the cascade failure. 0, 0| An abstract network 
is assumed for a power system in these models. The nodes in the network correspond to generators or consumers and 
the links correspond to transmission lines. Some removal rules of links in the network were assumed to simulate an 
elemental failure process. A cascade failure could be reproduced by numerical simulation of some simple models, and 
analytical theories were proposed. 

However, in actual alternating current (AC) systems, not only the network structure but also the electric potentials 
or phase variables are important to quantitatively describe a large-scale failure of power supply. For example, there 
is a phenomenon called voltage collapse, which occurs in power grids. In a voltage collapse, the electric potential at a 
consumer (a load) cannot be maintained at a stable stationary value and power supply to the load becomes impossible 
when the power demand is beyond a threshold, which can lead to a large-scale blackout. 

Taking the importance of phase variables into consideration, several authors proposed coupled oscillator models 
similar to the Kuramoto model as a model of power grids, and synchronization transition was studied in the models. [H- 
[Toj In synchronous generators, the generator turbine is driven by mechanical power, and the rotational motion is 
directly transformed into the sinusoidal oscillation of an alternating current. The phase variable of the alternating 
current for the ith generator is denoted as &{. The power produced at a generator is dissipated, accumulated, and 
transmitted along the electric transmission line. The output power is expressed using a sinusoidal function of the 
phase difference between the generators in the AC system. A simple energy balance equation leads to a set of phase 
equations: 



where = ftj — u>ot, ujq is the standard frequency, D is the damping constant, Wi is related to the input power 
supplied to the ith element, and F is the coupling constant. The phase variable <f>i denotes the difference between 
6i and the phase advance u)$t by the standard frequency. The authors performed some numerical simulations and 
obtained some theoretical results. 

When the damping effect is sufficiently large, the Kuramoto model is recovered as 



Phase transitions in globally coupled systems were studied in ref. 11 and ref. 12. Recently, synchronization in scale- 
free networks has been studied by several authors [l3l - fl5T ]. and a model with two subpopulations corresponding to 
generators and loads was also studied in ref. 16. In these models, the mutual synchronization is a main issue of research. 
As a result of the mutual synchronization, a macroscopic oscillation with a certain average frequency appears. The 
average frequency of the macroscopic oscillation is determined internally by self-organization in coupled oscillators, 
such as eqs. (1) and (2). 

However, it is important to maintain the standard frequency o>o at a constant value, such as 50 Hz or 60 Hz, in 
electric power systems. If the power demand exceeds a critical value, the standard frequency cannot be maintained, 
and the frequency begins to depart from the standard value. This is called loss of synchronization or step out. A step 
out seriously affects power grids. In this paper, we propose a new phase model of power grids, composed of sources 
(generators) and consumers (loads), that can simulate the occurrence of the step out phenomenon. We numerically 
study cascade failure in an artificially constructed networks: a square lattice and a scale-free network. The numerical 
results are compared with a mean-field approximation. 




(1) 




(2) 
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II. PHASE MODEL OF POWER GRIDS 



Details on electric power systems, including actual generators, can be found in textbooks on electrical engineering, 
e.g., ref. 17. In this paper, we study a very simplified model that includes generators and loads, although actual 
power systems are more complicated. For generators, we extend the phase model described by eq. (1). The equation 
of motion for is written as 

dt* ~ ° dt +Wl Fgz ' {d > 
where D is the damping constant, Wi is related to the input power, and P gi denotes the effective output power. Note 
that d(j>i /dt = d9i jdt — uj denotes the deviation of the frequency of the ith generator from the standard frequency 
uj . The output power P gi is expressed as 

Pgi = ^2/{EiEjYij sinOi - <f>j - a l} ) + E\ sin ay}, (4) 
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where Ei is the electric potential at the ith node, Yij = 1 / \J R\ ') + Wg Lf^ , and = tan -1 {Rij / (cooLij)} . Here 

denotes the resistance between the ith and jth nodes, and denotes the impedance between these two nodes. If 
Rij = 0, then = 0, and the sine of the phase difference <fii — (f>j multiplied by EiEj contributes to the effective power. 
The summation in eq. (4) is taken for all j nodes, including the generators and loads linked to the ith generator. In 
this paper, for simplicity, the output voltage Ei for generators is assumed to take a constant value Eq%- 

In contrast to the models studied in ref. 8 and ref. 9, the proposed model incorporates an additional equation to 
keep the effective frequency d<j)/dt at zero. We use a simple feedback control system expressed as 

dWi d(f> 

~df = ~ 7 di ' (5) 
This feedback control works as follows. If defy/dt is positive or the effective frequency is too fast, the input power is 
decreased. Conversely, if d<fi/dt is negative or the effective frequency is too slow, the input power is increased. If the 
feedback control works well, the effective frequency dfa/dt is maintained at zero. Under this condition, the entire 
power system, which includes the generators and loads, is synchronized at the standard frequency ujq, such as 50 Hz 
or 60 Hz. In actual power generators, the frequency is controlled by a local feedback using governors and a global 
regulation by the control center. Equation (5) corresponds to the local feedback by governors. Equation (5) implies 
that Wi(t) + j(j>i(t) =const. If the constant is set to zero, then Wi = —j<fii- Then, eq. (3) takes a form similar to that 
of the equation of motion of a damping harmonic oscillator with an output power term. 

d 2 <pi d<p t 

This is our model equation for generators. The form of the equation is similar to that of eq. (1); however, the 
dynamical behavior is rather different from that in the case of eq. (1). 

If the power demand expressed by P g i is too large, the input power Wi cannot respond to the excessive demand. If 
the generators become overloaded, they might break down. To avoid a breakdown, the generators might be stopped. 
In such an event, the overall synchronization would be lost; that is, a loss of synchronization or step out would occur, 
which might lead to a cascade failure. To express such a situation, we assume a simple rule that Ei is set to zero 
if Wi is beyond a critical value W c ,. If Ei(t) is set to zero, the function of the generator is lost. It is equivalent to 
effectively removing the generator node i from the power network system. 

For loads, on the other hand, we use different energy balance equations. That is, we assume that definite quantities 
of effective power P ei and reactive power Q ri are consumed at each load i. The complex power P ei +iQ ri is expressed as 
f^e 1 ^*/* , where Ii is the current supplied to the electric devices at the ith node, and * implies the complex conjugate. 
The current Ii is expressed as 

using the difference between the electric potentials and the admittance. Here, l/(Rij +iujLij) is denoted as — iYije laij . 
The powers P e i and Q ri are therefore expressed as 



Pei = ^ Yi 3^ EiE 3 sin (0J -<t>i+ a ij) - E i shlO^}, 
3 

Qri = ^ Yi 3'{ EiE 3 COS (4'3 ~<t>i + C^ij) - El cos a l} } . (7) 
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The summation is taken for all sites j, including the generators and loads linked to the ith node. 

Equations (7) are coupled algebraic equations. The phase <pi and electric potential Ei are determined for each i 
by solving eq. (7) for specified values of P e i and Q r i. We have numerically solved the coupled equations (7) with an 
iterative method. Generally, if P e i and Q r , are sufficiently small, both stable and unstable stationary solutions exist 
for Ei. However, if the effective power P e i is gradually increased as a control parameter, the magnitude of a stable 
solution Ei decreases and that of an unstable one increases. At a critical value of P e i, a transition occurs such that 
the stable and unstable solutions merge and disappear. There are no stationary solutions beyond the critical value. 
This phenomenon is called voltage collapse. If a voltage collapse occurs, the voltage Ei decreases rapidly, which can 
lead to an electric failure of the power supply. A voltage collapse occurs more easily if the reactive power Q ri is large. 
If a voltage collapse occurs, there is no stationary solution in the coupled equations (7). In our simple model, we have 
assumed another rule that Ei is set to zero when Ei(t) decreases to zero. This operation is equivalent to removing the 
load node i from the network of power grids, because Ei appears as a coefficient for the interaction terms in eq. (7). 

Thus, we have proposed a new simple phase model of power grids with two removal rules. The two removal rules 
simulate the situations of step out and the voltage collapse. Our model is similar to the models in ref. 6 and ref. 7 in 
that some removal rules are assumed. At the same time, our phase model has a form similar to that of eq. (1) studied 
in ref. 8 and ref. 9, although the dynamical behavior is rather different because feedback control is very important in 
our model. We perform numerical simulations of coupled equations (6) and (7) on a scale-free network for the sake of 
simplicity, although real power-grids are not regular square lattice and a scale-free network. Our model systems can 
be easily applied to other types of networks. Actually, we have performed numerical simulation on a square lattice 
and found cascade failure. The scale-free network is constructed using the preferential attachment method of Barabasi 
and Albert. New nodes are randomly assigned as being generators with probability p or loads with probability 1 — p. 

III. NUMERICAL SIMULATIONS ON A SCALE-FREE NETWORK 

In this paper, we show some numerical results for a uniform system where the parameters do not depend on z, i.e., 
P e i = P e , Q r i — Q r , Eci — Eq, W c i — W c , Yij — Y, and a t j — a. The total number of nodes is denoted by N. Some 
parameters in this paper are fixed, namely, D = 1,7 = l,Q r = 0.001, Y — 1, and Eg = 1, whereas N,W C , P e ,p and 
a are changed as control parameters. 

Figure 1 shows a numerical result of a scale-free network of N = 100 at p = 0.15, P e = 0.2, W c = 1.5, and a = 0. 
The initial conditions are (pi — 0,d4>i/dt — 0,Ei = 1, and W{ = 0. There are 15 generators and 85 loads in this 
network. Figure 1(a) shows a network at t = 2, where no failure occurs. The rhombi denote generator nodes and the 
plus signs denote load nodes. Figure 1(b) shows a network at t = 10, which is a stationary state. The input power 
Wi of one generator i = 6 of degree 14 goes beyond the critical value W c , and a loss of synchronization occurs; hence, 
the voltage Ei is set to zero. The generator node i = 6 is one of the hubs with large degree. The power demand to 
the other generators increases by the effective removal of this generator. As a result, a loss of synchronization occurs 
successively at two other generator nodes, namely, i = 1 of degree 24 and i — 9 of degree 13. Meanwhile, as a result of 
the shortage of the power supply, a large-scale voltage collapse occurs for the load nodes. When the voltage collapse 
occurs, the load nodes are effectively removed. Then, the successive loss of synchronization stops, because power 
demand to the other generators decreases. Finally, 12 generators and 5 loads survive in this numerical simulation, as 
shown in Fig. 1(b). As the numerical simulation suggests, the loss of synchronization tends to occur like a cascade 
failure; however, the voltage collapse tends to suppress the cascade failure because the power demand decreases when 
the load nodes are removed. 



IV. MEAN-FIELD APPROXIMATION 

Our system can be simplified with a mean-field approximation. The mean-field approximation is considered to be 
good if the degree of each node is sufficiently large. A similar mean-field analysis was performed for coupled phase 
oscillators in a scale-free network. [Til Il5| Note that D = 7 = 1 is assumed in the following equations. The phases 0i's 
of generators and loads at nodes of degree I are assumed to take the same value and are respectively denoted as (pGi 
and cf>Li. Similarly, the voltages Ei's at generator nodes and load nodes of degree I are assumed to take the same value 
and are respectively denoted as Eqi and En. The voltage Eqi takes a constant value Eq when the generator nodes 
do not exhibit step out, and it becomes zero when the generator nodes exhibit step out. On the other hand, E^s are 
determined by coupled equations (10) and (11) and become zero if the loads exhibit voltage collapse. The probability 
distribution of degree I is denoted as P(l). The distribution P{1) obeys a power law in a scale- free network, but our 
theory can be applied to any degree distribution P(l). We have used the distribution P(l) of the randomly constructed 
network used in our numerical simulation. \Y\ is denoted as Y. We introduce two kinds of order parameters, <jq and 
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FIG. 1. (a) Scale-free network at t = 2. Generators are denoted by rhombi and loads are denoted by +. (b) Network at t — 10. 
The generator nodes that exhibit step out and the load nodes that exhibit voltage collapse are removed; moreover, the links 
leaving the removed nodes are also removed. 



<jl, for generators and loads, respectively. 

J2P(l)lE Gl e^' 



a G exp(i<p G ) = 



a L exp(i(t) L ) = -£P(l)l ' 

They are average values of Ecie 1 ^ 01 and Eue 1 ^" weighted with degree Z. By using these order parameters, eq. (6) is 
approximated as 

'P&Gi , d(j>Gl , lv „ f . , , , , > _ . -. 

i>Gi TT" +plYE G i{a G sm((/) G - <pai + a) - E G i sin a} 



dt 2 ^" dt 

+ (1 - p)lYE Gl {a L sm(4> L ~ 4>Gl + a) - E G i sina}. (9) 
On the other hand, eq. (7) is approximated as 

P e = plYE u {a G sm(4> G - § Ll + a) - E u sina} 

+ (1 - p)lYE Ll {a L sinO L - <j> Ll + a) - E u sin a}, (10) 
Q r = plYE L i{cr G cos(4> G - (j>u + a) - E u cos a} 

+ (1 - p)lYE u {a L cos(c/) L - <f>u + a) - E u cos a}. 

(11) 

As stated already, the additional rule for step out is that E G \ is set to zero when Wi — —<f>Gi exceeds the critical value 
W c . The other additional rule for voltage collapse is that Eu is set to zero when there are no stationary solutions to 
eqs. (10) and (11). We have performed numerical simulation of coupled equations (8), (9), (10), and (11) to find a 
solution in the mean-field theory. 

Figure 2(a) compares the result for a scale-free network of N — 600 at a = 0,P e = 0.1, and p — 46/600 = 0.077 
obtained by direct numerical simulation with that obtained by mean-field approximation. The entire system is 
synchronized, and voltage collapse does not occur at these parameter values. The horizontal axis denotes degree Z, 
and the vertical axis denotes <j>i at the stationary state. The phases 4>i's for generators are denoted by +, and the 
0's for loads are denoted by rhombi. The two solid curves represent <p G i and 4>u obtained by numerical simulation of 
the mean-field equation. Rather good agreement is observed. As P e is increased, loss of synchronization occurs first, 
when Wi — —<\>i for a generator with a large degree Z exceeds the critical value W c = 1.5. The critical values P ec of 
P e for step out have been calculated as functions of p by direct numerical simulations of N = 600 and a = 0. The 
results are shown with rhombi in Fig. 2(b). The numerical result obtained using the mean-field equation is shown by 
a solid curve. Fairly good agreement is observed. It is natural that the critical value of the effective power increases 
with the number of the generators. 

The mean-field equation is still rather complicated. If we further assume that a G ~ 1, <tl *~ 1, a = 0, E G i — 1, Eu ~ 
X,Y = 1, and \4>gi ~ §u\ << 1 (which is approximately the case for small p), the critical value of step out can be 
approximately evaluated. The stationary solution of eq. (9) satisfies 



- 4>gi + pZ(0 G - 4>gi) + (i - p)1{4>l - 4>a) = o. 



(12) 
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FIG. 2. (a) Relation between the phase <f>i and the degree I in the scale-free network. The + signs and rhombi denote numerical 
results for generators and loads. Solid curves denote numerical results for the mean-field equation. The upper one corresponds 
to generators and the lower one corresponds to loads, (b) Critical values P ec as functions of the number ratio p of generators. 
The rhombi denote the results of the direct numerical simulation. The solid curve is the result obtained using the mean-field 
equation, and the dotted curve denotes eq. (17). (c) Number ratios vg and tl of removed nodes as functions of P e for p — 0.077. 



That is, 

plcpG + (1 - P)Hh 



l + l 



{P0 G + (1- J3 )0 L }(1-1//). (13) 



Similarly, eq. (10) yields 



9LZ 



pl<t> G + (1 - P)l<t>L ~ Pe 



(14) 



Substitution of eqs. (13) and (14) into eq. (8) yields 

exp{i(l -p)((/)G - 4>l)} 



-i{pij>G + (l-p)<pL}/l 



exp{ip(ct> L - «fc)} = ^ ^' ■ (15) 
By using the approximation exp(zx) ~ 1 + ix, the following approximate relations are obtained: 

4> G -<j>L=Pe/{p{l)), 
P<t>G + (1 ~P)H = -(1 ~P)Pe/p, 

<f>Gl = -(l-p)P e /p(l-l/l), 

Hi = -(1 - p)Pe/p ~ Pe/l, (16) 

where (I) is the average value of degree and is given by (/) = J P{l)ldl. As is seen from eq. (16), 4>gi decreases as l/l 
and 4>li increases as —l/l and 4>gi and <Pli approaches for I — > oo. The numerical results for the mean-field equation 
shown by the solid curves in Fig. 2(a) exhibit similar behavior that </>gi decreases as l/l and 4>u increases as — l/l for 
large I. The maximum value of Wi is evaluated as ~4>gi = (1 — p)Pe/p(l — l/l) ~ (1 ~p)Pe/p- Therefore, the critical 
value of P e is evaluated as 

P ec = f^- (17) 
1 ~P 

The dotted curve in Fig. 3(b) denotes eq. (17). The rough approximation seems to be sufficiently accurate for small 
P- 

We have numerically calculated tq and rz, also for the scale-free network of N = 600 . Figure 2(c) shows tq and 
Tl as functions of P e for p = 0.077. The solid curves denote results of the direct numerical simulation and the dashed 
curves denote results of the mean-field model. It is seen that the transitions are discontinuous, and the mean-field 
approximation is fairly good. 
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V. SUMMARY 



We have proposed a phase model of power grids. A simple feedback control was used to maintain the standard 
frequency in the model. Electric failure occurred via voltage collapse and step out (loss of synchronization). The 
two electric failures were incorporated as two removal rules for nodes of generators and loads. We performed direct 
numerical simulation and studied the mean-field model on a square lattice and in a scale-free network. We have 
found two kinds of transitions in a large network including many generators and loads. The step-out dominant failure 
expands as a cascade failure. The mean-field approximation was fairly suitable for understanding the complicated 
behavior of this system. In this paper, we used a regular square lattice and the scale-free network of Barabasi and 
Albert, and we assumed a uniform system for the sake of simplicity. Our model can be applied to general networks, 
and the mean-field theory might be applied to dense networks with large degrees. We have assumed that the critical 
value of power supply W c i and the power demand P e i are the same for all i's; however, in practice, these values are 
generally distributed widely. It is left to future research to study more realistic power grids with heterogeneous power 
supply and demand. 
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